%% 已知集热器坐标 D=156.0504 改变参数 W H
D=156.0504 ;
h=6;
% beta为同层相邻定日镜夹角
A=[];
for H=2:8
for W=H:8
%% 放置定日镜
dr=11+2*W/8;
R=100:dr:(350+D);
x=[];y=[];
theta_bond=zeros(1,length(R));
for i=1:length(R)
if R(i)<350-D % 定日最远边界
theta_bond(i)=-pi/2;
else
theta_bond(i)=asin((R(i)^2+D^2-350^2)/(2*R(i)*D));
end
end
beta=dr./R;% 满足相邻定日镜大于5m 约束6
for i=1:length(beta)
theta=( (-1)^i*beta(i)/4+theta_bond(i)+beta(i)/2):beta(i):(pi-theta_bond(i)-beta(i)/4);%更新theta
x=[x,R(i)*cos(theta)];
y=[y,R(i)*sin(theta)]; %放置定日镜
end
%% 计算本次的额定热功率 year_EF
location=[x',y'-D];
[year_EF,~]=Objfun2(H,W,h,D,location);
A=[A;W,H,year_EF];
end
end
%% H=6
A=[];
for H=6
for dr=(H+5):0.2:(H+7)
W=H;
R=100:dr:(350+D);
x=[];y=[];
theta_bond=zeros(1,length(R));
for i=1:length(R)
if R(i)<350-D % 定日最远边界
theta_bond(i)=-pi/2;
else
theta_bond(i)=asin((R(i)^2+D^2-350^2)/(2*R(i)*D));
end
end
beta=dr./R;% 满足相邻定日镜大于5m 约束6
for i=1:length(beta)
theta=( (-1)^i*beta(i)/4+theta_bond(i)+beta(i)/2):beta(i):(pi-theta_bond(i)-beta(i)/4);%更新theta
x=[x,R(i)*cos(theta)];y=[y,R(i)*sin(theta)]; %放置定日镜
end
location=[x',y'];
[year_EF,~]=Objfun2(H,W,h,D,location);
A=[A;W,dr,year_EF];
if size(A,1)>2 && A(end,3)-A(end-1,3)<0 && A(end-1,3)-A(end-2,3)<0
break
end
end
end
%% H=5
A=[];
for H=5
for dr=(H+5):0.2:(H+7)
W=H;
R=100:dr:(350+D);
x=[];y=[];
theta_bond=zeros(1,length(R));
for i=1:length(R)
if R(i)<350-D % 定日最远边界
theta_bond(i)=-pi/2;
else
theta_bond(i)=asin((R(i)^2+D^2-350^2)/(2*R(i)*D));
end
end
beta=dr./R;% 满足相邻定日镜大于5m 约束6
for i=1:length(beta)
theta=( (-1)^i*beta(i)/4+theta_bond(i)+beta(i)/2):beta(i):(pi-theta_bond(i)-beta(i)/4);%更新theta
x=[x,R(i)*cos(theta)];y=[y,R(i)*sin(theta)]; %放置定日镜
end
location=[x',y'];
[year_EF,~]=Objfun2(H,W,h,D,location);
A=[A;W,dr,year_EF];
if size(A,1)>2 && A(end,3)-A(end-1,3)<0 && A(end-1,3)-A(end-2,3)<0
break
end
end
end
%% H=4
A=[];
for H=4
for dr=(H+5):0.2:(H+7)
W=H;
R=100:dr:(350+D);
x=[];y=[];
theta_bond=zeros(1,length(R));
for i=1:length(R)
if R(i)<350-D % 定日最远边界
theta_bond(i)=-pi/2;
else
theta_bond(i)=asin((R(i)^2+D^2-350^2)/(2*R(i)*D));
end
end
beta=dr./R;% 满足相邻定日镜大于5m 约束6
for i=1:length(beta)
theta=( (-1)^i*beta(i)/4+theta_bond(i)+beta(i)/2):beta(i):(pi-theta_bond(i)-beta(i)/4);%更新theta
x=[x,R(i)*cos(theta)];y=[y,R(i)*sin(theta)]; %放置定日镜
end
location=[x',y'];
[year_EF,~]=Objfun2(H,W,h,D,location);
A=[A;W,dr,year_EF];
if size(A,1)>2 && A(end,3)-A(end-1,3)<0 && A(end-1,3)-A(end-2,3)<0
break
end
end
end
%% H=7
A=[];
for H=7
for dr=(H+5):0.2:(H+7)
W=H;
R=100:dr:(350+D);
x=[];y=[];
theta_bond=zeros(1,length(R));
for i=1:length(R)
if R(i)<350-D % 定日最远边界
theta_bond(i)=-pi/2;
else
theta_bond(i)=asin((R(i)^2+D^2-350^2)/(2*R(i)*D));
end
end
beta=dr./R;% 满足相邻定日镜大于5m 约束6
for i=1:length(beta)
theta=( (-1)^i*beta(i)/4+theta_bond(i)+beta(i)/2):beta(i):(pi-theta_bond(i)-beta(i)/4);%更新theta
x=[x,R(i)*cos(theta)];y=[y,R(i)*sin(theta)]; %放置定日镜
end
location=[x',y'];
[year_EF,~]=Objfun2(H,W,h,D,location);
A=[A;W,dr,year_EF];
if size(A,1)>2 && A(end,3)-A(end-1,3)<0 && A(end-1,3)-A(end-2,3)<0
break
end
end
end
%% H=8
A=[];
for H=8
for dr=(H+5):0.2:(H+7)
W=H;
R=100:dr:(350+D);
x=[];y=[];
theta_bond=zeros(1,length(R));
for i=1:length(R)
if R(i)<350-D % 定日最远边界
theta_bond(i)=-pi/2;
else
theta_bond(i)=asin((R(i)^2+D^2-350^2)/(2*R(i)*D));
end
end
beta=dr./R;% 满足相邻定日镜大于5m 约束6
for i=1:length(beta)
theta=( (-1)^i*beta(i)/4+theta_bond(i)+beta(i)/2):beta(i):(pi-theta_bond(i)-beta(i)/4);%更新theta
x=[x,R(i)*cos(theta)];y=[y,R(i)*sin(theta)]; %放置定日镜
end
location=[x',y'];
[year_EF,~]=Objfun2(H,W,h,D,location);
A=[A;W,dr,year_EF];
if size(A,1)>2 && A(end,3)-A(end-1,3)<0 && A(end-1,3)-A(end-2,3)<0
break;
end
end
end
e=0.05;
h=6;
mind=4;maxd=7;
leftd=(maxd-mind)/3;rightd=2*(maxd-mind)/3;
dddd=inf;
B=[];
while dddd>=e
leftd=mind+(maxd-mind)/3;rightd=mind+2*(maxd-mind)/3;
YEAR_EF=[];
for W=[leftd,rightd] % 左中右
H=W;
dr=W+5;
R=100:dr:(350+D);
x=[];y=[];
theta_bond=zeros(1,length(R));
for i=1:length(R)
if R(i)<350-D % 定日最远边界
theta_bond(i)=-pi/2;
else
theta_bond(i)=asin((R(i)^2+D^2-350^2)/(2*R(i)*D));
end
end
beta=dr./R;% 满足相邻定日镜大于5m 约束6
for i=1:length(beta)
theta=( (-1)^i*beta(i)/4+theta_bond(i)+beta(i)/2):beta(i):(pi-theta_bond(i)-beta(i)/4);%更新theta
x=[x,R(i)*cos(theta)];y=[y,R(i)*sin(theta)]; %放置定日镜
end
% 计算本次的额定热功率 year_EF
location=[x',y'-D];
h=6;
[aaa,~]=Objfun2(H,W,h,D,location);
YEAR_EF=[YEAR_EF,aaa];
end
if YEAR_EF(1)<YEAR_EF(2)
B=[B,mind,leftd,rightd,maxd,YEAR_EF];
mind=leftd;
dddd=abs(leftd-rightd)
elseif YEAR_EF(1)>YEAR_EF(2)
B=[B,mind,leftd,rightd,maxd,YEAR_EF];
maxd=rightd;
dddd=abs(leftd-rightd)
end
end
save('D','D')
save('W','W')



%3
clear;
% load D;load W;
W=5.4670;
D=237.0240;
dr=W+5;
h0=6;
H=W;
R=100:dr:(350+D);
x=[];y=[];num=[];
theta_bond=zeros(1,length(R));
for i=1:length(R)
if R(i)<350-D % 定日最远边界
theta_bond(i)=-pi/2;
else
theta_bond(i)=asin((R(i)^2+D^2-350^2)/(2*R(i)*D));
end
end
beta=dr./R;% 满足相邻定日镜大于5m 约束6
for i=1:length(beta)
theta=( (-1)^i*beta(i)/4+theta_bond(i)+beta(i)/2):beta(i):(pi-theta_bond(i)-beta(i)/4);%更新theta
x=[x,R(i)*cos(theta)];y=[y,R(i)*sin(theta)]; %放置定日镜
num=[num,length(x)];
end
location=[x',y'-D];
n=length(beta);
step=(6-W/2)/n;
h_tj0=W/2+step/2:step:6;
maxval=0;
for time=1:1e6
h_tj=h_tj0+step*rand(1,length(beta))-step/2; %随机取点
% h_tj=sort(h_tj);
h(1:num(i))=h_tj(1); % h 赋值过程
for i=2:length(beta)
h(num(i-1)+1:num(i))=h_tj(i);
end
[year_Ef,year_Ef_per_area]=Objfun3(H,W,h,D,location)
if year_Ef>maxval
maxval=year_Ef;
maxh=h;
maxh_tj=h_tj;
end
end



